title: “presentation_happign_fif” author: “Paul CARTERON” format: html
df-print: paged
self-contained: trueeditor: visual
library(happign)## Please make sure you have an internet connection.## Use happign::get_last_news() to display latest geoservice news.library(sf)## Linking to GEOS 3.11.2, GDAL 3.8.2, PROJ 9.3.1; sf_use_s2() is TRUElibrary(terra)## terra 1.7.78library(tmap);tmap_mode("view")##
## Attaching package: 'tmap'## The following object is masked from 'package:datasets':
##
## rivers## tmap mode set to 'view'Pour obtenir une ressource à partir des flux WFS ou WMS de l’IGN il faut :
get_wms_raster() ou
get_wfs()sf::st_read() ou de facon interactive avec
mapedit::drawFeatures()get_layers_metadata()shp <- read_sf(system.file("extdata/penmarch.shp", package = "happign"))
tm_shape(shp)+
tm_polygons("red")get_apikeys()## [1] "administratif" "adresse" "agriculture" "altimetrie"
## [5] "cartes" "cartovecto" "clc" "economie"
## [9] "enr" "environnement" "geodesie" "lambert93"
## [13] "ocsge" "ortho" "orthohisto" "parcellaire"
## [17] "satellite" "sol" "topographie" "transports"apikey <- get_apikeys()[1]layers <- get_layers_metadata("wfs", "administratif")
layer <- layers[1,1]arrondissement <- get_wfs(x = shp,
layer = "ADMINEXPRESS-COG-CARTO.LATEST:arrondissement")## Features downloaded : 1tm_shape(arrondissement)+
tm_polygons()+
tm_shape(shp)+
tm_polygons("red")apikey <- "parcellaire"
layer <- "CADASTRALPARCELS.PARCELLAIRE_EXPRESS:parcelle"
get_layers_metadata("wfs", "parcellaire")## Name
## 1 BDPARCELLAIRE-VECTEUR_WLD_BDD_WGS84G:arrondissement
## 2 BDPARCELLAIRE-VECTEUR_WLD_BDD_WGS84G:batiment
## 3 BDPARCELLAIRE-VECTEUR_WLD_BDD_WGS84G:commune
## 4 BDPARCELLAIRE-VECTEUR_WLD_BDD_WGS84G:divcad
## 5 BDPARCELLAIRE-VECTEUR_WLD_BDD_WGS84G:localisant
## 6 BDPARCELLAIRE-VECTEUR_WLD_BDD_WGS84G:parcelle
## 7 AOC-VITICOLES:aire_parcellaire
## 8 CADASTRALPARCELS.PARCELLAIRE_EXPRESS:arrondissement
## 9 CADASTRALPARCELS.PARCELLAIRE_EXPRESS:borne_limite_propriete
## 10 CADASTRALPARCELS.PARCELLAIRE_EXPRESS:batiment
## 11 CADASTRALPARCELS.PARCELLAIRE_EXPRESS:commune
## 12 CADASTRALPARCELS.PARCELLAIRE_EXPRESS:feuille
## 13 CADASTRALPARCELS.PARCELLAIRE_EXPRESS:borne_parcelle
## 14 CADASTRALPARCELS.PARCELLAIRE_EXPRESS:localisant
## 15 CADASTRALPARCELS.PARCELLAIRE_EXPRESS:parcelle
## 16 CADASTRALPARCELS.PARCELLAIRE_EXPRESS:subdivision_fiscale
## Title
## 1 BDParcellaire vecteur (WFS)
## 2 BDParcellaire vecteur (WFS)
## 3 BDParcellaire vecteur (WFS)
## 4 BDParcellaire vecteur (WFS)
## 5 BDParcellaire vecteur (WFS)
## 6 BDParcellaire vecteur (WFS)
## 7 Délimitation parcellaire AOC viticole
## 8 PCI Vecteur : Arrondissements
## 9 PCI Vecteur : Bornes de limite de propriété
## 10 PCI Vecteur : Bâtiments
## 11 PCI Vecteur : Communes
## 12 PCI Vecteur : Feuilles
## 13 PCI Vecteur : Lien bornes parcelles
## 14 PCI Vecteur : Localisants
## 15 PCI Vecteur : Parcelles
## 16 PCI Vecteur : Subdivisions fiscales
## Abstract
## 1 BDParcellaire vecteur.
## 2 BDParcellaire vecteur.
## 3 BDParcellaire vecteur.
## 4 BDParcellaire vecteur.
## 5 BDParcellaire vecteur.
## 6 BDParcellaire vecteur.
## 7 Selon les règlements européens n°110/2008 du 15 janvier 2008, n°1151/2012 du 21 novembre 2012 et n°1308/2013 de l’OCM , le cahier des charges des appellations définit l’aire géographique. Le terme de « délimitation parcellaire » désigne une aire qui repose sur les limites administratives du cadastre (les parcelles) et dont le maillage suffisamment fin permet de tenir compte de variations très localisées des éléments du milieu physique. Elle correspond à l'aire de production de la matière première. Elle est incluse dans l'aire géographique. Données produites par l'Institut National des Origines et la Qualité (INAO)
## 8 Arrondissements avril 2023
## 9 Bornes de limite de propriété avril 2023
## 10 Bâtiments avril 2023
## 11 Communes avril 2023
## 12 Feuilles avril 2023
## 13 Lien bornes parcelles avril 2023
## 14 Localisants avril 2023
## 15 Parcelles avril 2023
## 16 Subdivisions fiscales avril 2023# bbox
parcel_bbox <- get_wfs(shp, layer, spatial_filter = "bbox")## Features downloaded : 1000...1303tm_shape(shp)+
tm_borders("black", lwd = 3)+
tm_shape(st_as_sfc(st_bbox(shp)), name = "bbox", is.master = T)+
tm_borders("red", lwd = 2)+
tm_shape(parcel_bbox)+
tm_borders("grey30")# intersects
parcel_intersects <- get_wfs(shp, layer, spatial_filter = "intersects")## Features downloaded : 1000...1110tm_shape(shp)+
tm_borders("red", lwd = 2)+
tm_shape(parcel_intersects)+
tm_borders("grey30")# within
parcel_within <- get_wfs(shp, layer, spatial_filter = "within")## Features downloaded : 983tm_shape(shp)+
tm_borders("red", lwd = 3)+
tm_shape(parcel_within)+
tm_borders("grey30")# dwithin
parcel_dwithin <- get_wfs(shp, layer,
spatial_filter = c("dwithin", 200, "meters"))## Features downloaded : 1000...1758tm_shape(shp)+
tm_borders("red", lwd = 3)+
tm_shape(st_buffer(shp, 200))+
tm_borders("red", lwd = 1)+
tm_shape(parcel_dwithin)+
tm_borders("grey30")Les filtres ECQLs sont une version simplifiée des filtres SQL. Il
permettent de filtrer sur l’ensemble de la source donnée de façon
précise. Pour utiliser un filtre ECQL il est nécessaire de connaitre le
nom des attributs de la table. Pour cela on utilise la fonction
get_wfs_attributes.
layer <- "LIMITES_ADMINISTRATIVES_EXPRESS.LATEST:commune"
attributs <- get_wfs_attributes(layer)
print(attributs)## [1] "id" "nom" "nom_m" "insee_com" "statut"
## [6] "population" "insee_can" "insee_arr" "insee_dep" "insee_reg"
## [11] "siren_epci"Une fois que les noms d’attributs sont connus, il suffit de construire le filtre ECQL. Dans cet exemple, on filtre toutes les communes commençant par F, terminant par F et contenant un I.
** Remarque : pour pouvoir requêter sur l’ensemble de la table
(france entière), il ne faut pas fournir une shape à
get_wfs i.e. x = NULL**
ecql_filter <- "nom_m LIKE 'F%I%F'"
fif_commune <- get_wfs(x = NULL,
layer,
ecql_filter = ecql_filter)## Features downloaded : 4tm_shape(fif_commune)+
tm_borders("red")+
tm_text("nom_m", col = "red")Il est également possible de combiner plusieurs filtres ECQLs et des filtres spatiaux. Dans cet exemple, toutes les communes précédente avec une population de plus de 1000 habitants.
ecql_filter <- "nom_m LIKE 'F%I%F' AND population > 1000"
fif_commune_plus_de_1000 <- get_wfs(x = NULL,
layer,
ecql_filter = ecql_filter)## Features downloaded : 1tm_shape(fif_commune_plus_de_1000)+
tm_borders("red")+
tm_text("nom_m", col = "red")res <- get_wfs(shp, interactive = T)
# qtm : quick tmap : permet d'afficher sans connaitre le type de donnée en amont
tm_shape(shp)+
tm_borders()+
qtm(st_make_valid(res))Pour récupérer un raster, le principe de base est le même que pour
get_wfs(), il faut uniquement préciser la résolution en
plus.
Dans cet exemple, le MNT de Penmarch est téléchargé à une résolution de 25m.
apikey <- "altimetrie"
layer <- "ELEVATION.ELEVATIONGRIDCOVERAGE" #get_layers_metadata(apikey, "wms")
penmarch <- get_wfs(shp, "LIMITES_ADMINISTRATIVES_EXPRESS.LATEST:commune")## Features downloaded : 1mnt <- get_wms_raster(x = penmarch,
layer = layer,
res = 25,
rgb = FALSE)## 0...10...20...30...40...50...60...70...80...90...100 - done.## Raster is saved at : C:\Users\PaulCarteron\AppData\Local\Temp\Rtmp2jjH9w\file5ed03a406ea9.tiftm_shape(mnt)+
tm_raster()+
tm_shape(penmarch)+
tm_borders("black", lwd = 2)Il est possible de demander un système de coordonnée ce qui évite de dévoir reprojeter le raster après téléchargement (chronophage). Cependant, il faut faire attention à préciser la résolution dans le système de coordonnée souhaitée. Par exemple pour le WGS84 (EPSG:4326), l’unité de mesure est le degré.
mnt <- get_wms_raster(x = penmarch,
layer = "ELEVATION.ELEVATIONGRIDCOVERAGE",
res = 0.0005,
crs = 4326,
rgb = FALSE)## 0...10...20...30...40...50...60...70...80...90...100 - done.## Raster is saved at : C:\Users\PaulCarteron\AppData\Local\Temp\Rtmp2jjH9w\file5ed047d66af7.tiftm_shape(mnt)+
tm_raster()+
tm_shape(penmarch)+
tm_borders("black", lwd = 2)Les API Carto sont mise à disposition par l’IGN pour faciliter le
remplissage de formulaire. Elles ont l’avantages de ne pas
nécessairement avoir besoin d’une shape en entrée. happign permet de se
connecter à ces APIs à l’aide des fonctions get_apicarto_*.
Toutes ces fonctions sont vectorisées, ainsi chaque paramètre peut
prendre plusieurs valeurs.
Données récupérables à parti du code insee, code département ou d’une shape.
com <- get_apicarto_cadastre(c("29158", "29165"), type = "commune")## Features downloaded : 1## Features downloaded : 1# com <- get_apicarto_cadastre("29158", type = "parcelle", section = "AB", numero = "0001")
tm_shape(com)+
tm_polygons()Permet de récupérer le code insee et nom d’une commune à partir du code postal.
code_post <- get_apicarto_codes_postaux(c("29760", "29260"))API connectée au géoportail de l’urbanisme permettant de récupérer divers documents d’urbanisme, prescriptions et servitudes d’utilités publiques (Cône de vue, Alignement d’arbre, Canalisations de transport de gaz, Espace boisé classé, Carrières, Périmètres monuments historiques, …)
Carte des Espaces Boisés Classés et des alignements d’arbres :
prescriptions <- get_apicarto_gpu("DU_93014",
ressource = c("prescription-surf", "prescription-lin"))## Features downloaded : 159ebc <- prescriptions[prescriptions$libelle == "EBC", ]
align_arbre <- prescriptions[prescriptions$libelle == "Alignement d'arbre", ]
com <- get_apicarto_cadastre("93014", "commune")## Features downloaded : 1tm_shape(com)+
tm_borders(lwd = 2)+
tm_shape(ebc)+
tm_fill("firebrick")+
tm_shape(align_arbre)+
tm_lines("purple")Remarque : Bien souvent, les APIcarto ne prennent pas ne charge les géometries trop complexes. La paramètre dTolerance permet de simplifier temporairement la géométrie.
penmarch <- get_apicarto_cadastre("29158", type = "commune")## Features downloaded : 1rpg <- get_apicarto_rpg(penmarch, 2020, dTolerance = 10) |>
st_make_valid()## Features downloaded : 141tm_shape(penmarch)+
tm_borders(lwd = 2)+
tm_shape(rpg)+
tm_polygons("code_cultu")Indisponible pour le moment !
Dans cet exemple, le MNT et le MNS permettent de calculer le MNH. Les courbes de niveau sopnt également rajoutées.
shp <- get_apicarto_cadastre("35211", "commune")## Features downloaded : 1layers <- get_layers_metadata("wms-r", "altimetrie")
mnt_layer <- layers[3,1]
mns_layer <- layers[4,1]
mnt <- get_wms_raster(shp, mnt_layer, res = 10, crs = 2154, rgb = FALSE, verbose = FALSE)## 0...10...20...30...40...50...60...70...80...90...100 - done.mns <- get_wms_raster(shp, mns_layer, res = 10, crs = 2154, rgb = FALSE, verbose = FALSE)## 0...10...20...30...40...50...60...70...80...90...100 - done.# Calculate digital height model i.e. tree height
mnh <- mns - mnt
mnh[mnh < 0] <- NA # Remove negative value
level_curve <- get_wfs(shp, "ELEVATION.CONTOUR.LINE:courbe",
spatial_filter = "intersects")## Features downloaded : 151level_curve_intersect <- st_intersection(level_curve, shp)## Warning: attribute variables are assumed to be spatially constant throughout
## all geometriestm_shape(mnh) +
tm_raster(col.scale = tm_scale_continuous(value.na = "grey", values = "-spectral"),
col.legend = tm_legend("Height")) +
tm_shape(level_curve_intersect)+
tm_lines(col = "black")+
tm_shape(shp)+
tm_borders(lwd = 2, col = "red")